Multi-scale image normalization and enhancement

ABSTRACT

An image may be processed to normalize and/or remove noise from the image. The processing of the image may involve decomposition of the image into multiple components and subsequent gray scale registration across multiple scales.

BACKGROUND

Radiographic imaging for medical purposes is well known in the art.Radiographic images of the chest, for example, provide importantdiagnostic information for detecting and treating a large number ofmedical conditions involving the lungs, bony structures in the chest,the upper abdominal organs, the vascular structures of the lungs, andthe disc spaces of the mid-thoracic spine.

Because of the great advantages provided by digital images, radiographsare increasingly stored and manipulated in digital form. Digitalradiographs may be created either by direct capture of the originalimage in digital form, or by conversion of an image acquired by an“analog” system to digital form. Digital images simplify record keeping,such as in matching radiographs to the correct patient, and allow formore efficient storage and distribution. Digital images also allow fordigital correction and enhancement of radiographs, and for applicationof computer-aided diagnostics and treatment.

Radiologists are highly skilled in interpreting radiographic images, butlimitations of radiographic systems, and variability between systems,can hamper proper interpretation. Sources of variability related to theacquisition of radiographic images may include the spatial sampling ofthe images; the gray scale resolution (or “bit depth”); the ModulationTransfer Function (“MTF”) of the system; image contrast; and noise.

The sampling function of an image can generally be expressed as thenumber of pixels in a unit length. Generally, sampling is performed ator near the Nyquist rate to avoid aliasing. For example, ahighly-detailed chest radiogram may have 5,000 pixels per inch, for aminimum discernable feature size of 200 microns. A uniform spatialresolution between images can be important in automated systems, such aswhen software is used to identify or analyze features having specificspatial characteristics in a radiogram.

Bit depth is the number of data bits used to store the brightness valueof each pixel of an image. Different radiographic systems may produceradiograms with different bit depths. For example, bit depths commonlyrange from 10 to 12 bits. Bit depth is important not just in respect tothe quality of the original image, but becomes a limiting factor whendigitally manipulating images, such as when processing the images toaccentuate particular features or in computer-aided diagnosis.Insufficient bit depth can result in degraded processed images, imagingartifacts, and unreliable diagnostic results.

Modulation Transfer Function (MTF) is the spatial frequency response ofan imaging system, or of an imaging component. High spatial frequenciescorrespond to fine image detail, while low spatial frequenciescorrespond to larger structures. The contrast produced on a radiographicimage by features of different sizes may differ due to the system MTF.Typically, the contrast of features at a high spatial frequency can bereduced relative to the contrast of features at a low spatial frequencydue to the limited resolving power of the imaging instrument. Because ofthe reduction in amplitude variation of smaller features due to MTF, thevisibility of smaller features in a radiograph may be masked byoverlying larger structures in the image.

Contrast involves the brightness differences between neighboring pixelsin an image. Contrast concerns not just the absolute difference betweenthe brightest and darkest pixels, but also the brightness distributionof the intermediate pixels. For example, the distribution of brightnessvalues may be skewed towards the bright or dark end of the distributionrange, making it difficult to discern features having similarbrightness. For both the human observer and for automated systems, it isbeneficial that different radiographic images have substantially similarcontrast to ensure consistent interpretation during reading orprocessing.

One technique used to correct for differences in gray scale appearanceis known as histogram matching. A histogram is essentially a bar graphrepresentation of the distribution of pixel values in an image, in whichthe heights of the bars are proportional to the number of pixels in theimage having that pixel value. As is known in the art, histogrammatching is a pixel mapping derived from an input cumulative densityfunction (CDF) and a target CDF. As CDFs are monotonic and lie in the0-1 range, histogram matching is a simple matter of alignment, ormatching, the two CDFs. Histogram matching is generally used in areaswhere it is of interest to be able to directly compare pixels values ofsimilar scenes; it is a global technique in that it uses pixel valuesand not spatial information in anyway. It does not, for example, addressthe problem localized contrast differences in images of the same scene.

Noise is a universal limitation of all measurement systems, includingradiographic systems. In radiographic systems, it is generally necessaryto limit the cumulative exposure of a subject to x-rays; a tradeoff forshort exposure times is an increase in noise on the resulting image.Noise on radiographic images tends to primarily manifest themselves athigher spatial frequencies. At the highest spatial frequencies, noisemay predominate, limiting the ability to discern detail in an image.

Both to improve the uniformity of radiographic images for interpretationby radiologists and other professionals, and to provide a goodfoundation for subsequent digital processing and analysis of the images,a principled system for normalization of images may be desirable. Thenormalization process may account for differences between images atdifference spatial frequencies and may address the problems of contrastand noise.

SUMMARY

Embodiments of the invention may include methods that convertradiographic images, such as chest radiograms, to a uniform size, orpixel spacing, and bit depth; decompose the images into spatialcomponents; separately adjust the gray-scale distribution of the spatialcomponents; reduce image noise; and produce a re-combined enhancedoutput image. The gray-scale distribution of the spatial components maybe based on a statistical ensemble, derived from the source of theradiographic images, an analytic model, or derived on a per image basisto maximize or minimize an objective function.

Various embodiments of the invention may be in the forms of methods,apparatus, systems, and/or computer-readable media containingprocessor-executable instructions to execute methods. It is furthernoted that it is anticipated that such methods may be performed by anautomated processing device, for example, but not limited to, an imageprocessing device, by a general purpose processor or computer, etc.

BRIEF DESCRIPTION OF THE DRAWINGS

FIG. 1 is a flow diagram providing an overview of embodiments of theinvention;

FIG. 2 further illustrates the resizing of an image, with FIG. 2( a)showing a portion of an original image at a first resolution, and FIG.2( b) showing the a portion of the resized image at a second resolution;

FIG. 3 is a flow diagram further illustrating the bit-depth adjustmentof an image according to an embodiment of the invention;

FIG. 4 schematically illustrates bit depth adjustment at a histogramlevel, with FIG. 4( a) illustrating a histogram of an original image,and FIG. 4( b) illustrating the histogram of the bit-depth adjustedimage;

FIG. 5 is a flow diagram illustrating multi-scale decomposition of animage according to an embodiment of the invention;

FIG. 6 shows for illustrative purposes a representative portion of achest radiograph;

FIG. 7 further illustrates the multi-scale decomposition of an image,with FIG. 7( a) illustrating the highest spatial frequencies of thedecomposed image; FIG. 7( b) illustrating the next lower spatialfrequency components; FIG. 7( c) illustrating the third spatialfrequency components; FIG. 7( d) illustrating the lower spatialfrequency components from the decomposition; and FIG. 7( e) illustratingthe “residual” image after decomposition; and

FIG. 8 illustrates a conceptual block diagram of a system in which allor a part of various embodiments of the invention may be implemented.

DESCRIPTION OF EMBODIMENTS OF THE INVENTION

FIG. 1 is a flow diagram providing an overview of embodiments of theinvention. Embodiments of the invention may typically begin with animage 102, such as a chest x-ray. The various systems used to “capture”x-rays may produce images having different resolutions (pixels per unitarea) and different bit depths (the number of data bits used to storeeach pixel). The method may thus begin with resizing 104 the raw image,if necessary, to a standard resolution for subsequent processing. Anexemplary embodiment of the invention, for example, utilizes an imageresolution of 5,000 pixels per linear inch, for a resolution of 200microns.

To provide a standard bit depth for subsequent processing, the image bitdepth may then be adjusted 106. The limited bit depth of some capturesystems may not be adequate for the subsequent processing steps in thenormalization process, and might result in a poor image quality orprocessing artifacts. An exemplary embodiment of the invention uses abit depth of 10 bits per pixel, allowing for 1024 intensity levels, asfurther explained with respect to FIGS. 3 and 4, below.

Although shown in FIG. 1 as occurring after the resizing of the image,in other embodiments bit depth adjustment may be performed first.

After the image has been modified to a standard resolution and bitdepth, a multi-scale (or multi-resolution) decomposition 108 may beperformed to separate the image into sub-images, which may containdifferent frequency components from the original image. There are avariety of methods that may be employed for multi-scale decomposition,including Steerable or Laplacian pyramids, wavelets, curvelets, and soon. An exemplary embodiment of the invention, to which the invention isnot limited, utilizes a non-decimated wavelet decomposition with aB3-spline as the scaling function.

Typically, multi-scale decompositions may result in sub-images coveringa frequency range of 2:1, with the sub-image including the highestfrequency components covering the upper half of the frequencies; andwith each subsequent sub-image covering the upper half of the remainingbandwidth. A “residual” image may also be generated, which covers thelower frequencies not otherwise included in the other sub-images.

Multi-scale decompositions may be “decimating”, in that each subsequentlevel of the decomposition may result in an image smaller in size thanan image of the previous level (for example, each subsequent level mayhave half the number of pixels, both horizontally and vertically).Decimation may allow for more efficient algorithms and therefore fasterexecution. Since it may be useful, however, to have all the sub-imagesof the same size and resolution in the gray-scale registration stepdiscussed below, exemplary embodiments of the invention may utilizenon-decimating decomposition.

In an exemplary embodiment, to which the invention is not limited, theimage is decomposed into a multi-scale representation using a redundant(non-decimating) wavelet transform, resulting in 7 levels, with thehighest level corresponding to spatial frequencies up to the 200 micronpixel spacing and the lowest level corresponding to approximately 13 cm.The decomposition in this exemplary embodiment may be non-decimating, inorder achieve translation invariance and potential artifacts associatedwith reconstruction. The exemplary embodiment utilizes a B3-spline asthe scaling function (smoothing filter), which may have band-passcharacteristics similar to that of a Gaussian kernel; the decompositionmay, therefore, be similar to a Laplacian pyramid. The decomposition ofthe image is discussed further with respect to FIGS. 5, 6 and 7, below.

After decomposition, the resultant sub-images may be subjected togray-scale registration 110 and noise removal and enhancement 112. Thesub-images from the decomposition of the image may each capturestructure for a particular scale. In order to account for sensorvariation associated with MTF and contrast, each sub-image may be mappedto a target distribution using histogram specification, as known in theart. In effect, embodiments of the invention may combine the localnature of the decomposition with the global matching properties ofhistogram matching.

In an exemplary embodiment, the target distributions may be defined byCumulative Density Functions (CDFs) that may be formed by averagingmultiple images from a set of radiographic images. In an exemplaryembodiment, all data may be taken from the lung-field, which maytypically include both air-filled and opaque portions. It is envisionedthat this averaging process may serve to reduce the patient specificityof images while reinforcing the persistent sensor characteristics, andto thereby minimize irrelevant “noise” associated with a particularacquisition device and/or patient.

Continuing with the description of embodiments of the invention, FIG. 2illustrates an example of the resizing of the raw image in additionaldetail. As seen in FIG. 2( a), a “raw” image 202 may be composed ofpixels 204 having a first resolution. To convert the raw image to anoutput image 212 having pixels 214 of a standardized second resolution,one of many known resizing algorithms may be used, including “nearestneighbor”, “bilinear”, and “bicubic”. An exemplary embodiment of theinvention, to which the invention is not limited, uses bilinearinterpolation, which may be used to linearly interpolate the value of apixel from the values of the nearest 4 original pixels.

FIG. 3 is a flow diagram further illustrating the bit-depth adjustmentof the original raw image 302 according to an embodiment of theinvention, and FIG. 4 schematically illustrates at a histogram level howbit depth adjustment may be accomplished. As seen at 402, the raw inputimage may not fully use the dynamic range available in its nativeformat, in that the ends of the corresponding histogram may not extendfully to end of the “dark” and “light” ranges. As better seen at 404,the pixel grayscale values may range from a first arbitrary value “N” toanother arbitrary value “M”, rather than from the lowest possiblegrayscale value (i.e., zero) to the highest possible grayscale value.

In an exemplary embodiment of the invention, bit depth adjustment mayfirst map 304 the minimum grayscale value of the original image (denoted“N” at 404) to, for example, zero in the resized image, and may then map306 the maximum grayscale value of the original image (denoted “M” at404) to, for example, 1023, the highest value possible with a bit depthof 10 (an example to which the invention is not limited). The remainingpixel values of the original image, in this example, may then beuniformly distributed 308 between “1” and “1022” in the bit-depthadjusted image 310, which may result in a distribution such as shown at414 (in increasing the bit depth, it may be noted that not all grayscalevalues in the adjusted image may be utilized).

The mapping of the minimum grayscale value in the original image to “0”and the maximum grayscale value to “1023”, while limiting the remaininggrayscale values to the range of 1 to 1022, may serve to preserve theminimum and maximum, which can later help with artifact detection.Again, it is noted that the mapping from 0 to 1023 is an example towhich the invention is not limited.

FIGS. 5, 6, and 7 further illustrate the multi-scale decomposition of aresized and bit-depth adjusted image according to an exemplaryembodiment of the invention. FIG. 5 is a flow diagram illustrating anembodiment of image decomposition. A non-decimating 7-level multi-scalewavelet decomposition 504 may be performed on the resized and bit depthadjusted image 502 (it is noted that such the invention is not limitedto this particular decomposition, however), which may result insub-images containing an upper band of frequency components 510,intermediate bands (not shown), and a lower band of frequency components520. A residual image 530, including those spatial frequencies below thelower band of the decomposition, may also be created 506.

FIG. 6 shows, for illustrative purposes only, a portion of ahypothetical chest radiogram, while FIG. 7( a) shows the hypotheticalradiogram depicted as an array of pixel values (it should be noted thatthe resolution does not correspond to that of a true radiogram, and thatthe figures serve only as a tutorial device). The bright areas in FIG. 6can be seen to correspond to high pixel values in FIG. 7, and the darkerareas of FIG. 6 to lower pixel values.

FIG. 7( b) illustrates the highest frequency components from an imagedecomposition of the hypothetical radiogram. It may be observed thatnoise contributes largely to the pixel values. In a true radiographicimage, it has been observed that there may often be little of interestin the highest level of the decomposition, which may mostly comprisenoise.

FIG. 7( c) illustrates the next level of a decomposition of thehypothetical radiogram. It may be observed that more actual structurefrom the image is apparent. FIG. 7( d) shows the lowest level of adecomposition of the hypothetical radiogram (levels intermediate between7(c) and 7(d) may also exist but have been omitted); it may be observedthat larger structures, having lower spatial frequencies, may bevisible. FIG. 7( e) shows the “residual” image from a decomposition ofthe hypothetical radiogram, which is what remains of the original imageafter the decomposition levels are removed.

The multi-scale decomposition, gray scale registration, and noiseremoval processes may be carried out jointly, according to variousembodiments of the invention. Each multi-scale detail may besuccessively generated and processed for noise removal, and gray scaleregistration may be accomplished. In one exemplary embodiment of theinvention, the only noise removal that occurs may be to leave the firstmulti-scale detail out of the reconstruction. For chest radiographs,this detail may often contain very little information and may be almostentirely noise.

In order to register the gray scale values, each multi-scale detail maybe subjected to a model matching process. The model may be derivedempirically or analytically and may be based on a-priori knowledge of atarget distribution or derived “on-the-fly” to achieve a desired resultsuch as maximum signal-to-noise for target scale(s). This may be used tomap the multi-scale details to a target distribution, may suppress andenhance the overall content at each scale, and may be used to accountfor variations in contrast, sharpness, and/or brightness, which may thusallow the method to operate across a wide variety of acquisitionsettings. By successively adding these registered details, a normalizedimage may be formed. Such a capability may allow known sensor artifactsto be suppressed while retaining the target signal.

The residual image of the multi-scale transform may be kept separatefrom the normalized part. The normalized image may represent thestructural content of the image, while the coarse image may representthe low-frequency content that is typically patient and/or dose specificand may carry very little information.

One may, thereby, obtain two images: one image, the normalized image,may correspond to the reconstructed multi-scale details that have beennormalized; the other image may correspond to a low-pass residual thatmay only contain gross global differences in image (this image, whilenot necessarily being included in all subsequent processes, may be addedback at the end to preserve the relative appearance of different areas;it is also noted that this component may be dynamically weighted tothereby provide different degrees of tissue equalization).

The degree of visual enhancement may be adapted based on an explicitsegmentation of region(s) or adapted based on an implicit estimate ofregion's density. In chest X-rays it is well known that the latitude ofthe image may be quite large due to the wide range of absorptionproperties that may exist within the chest. For example, the heart anddiaphragm may generally have a high absorption coefficient and may,therefore, require far greater enhancement than the air-filled region ofthe lung.

In the case of dynamic (implicit) adaptation to density, according tosome embodiments of the invention, the following may be performed:Enhanced image=G*C+K1*R,  [eq 1]where C may represent a registered/enhanced contrast component from thenormalization; R may represent the coarse residual; and K1 is a scalardesign parameter for dynamic range reduction; furthermore, G may begiven by:G=1+(K2)*{circumflex over (R)}^(K3)−(K4)*{circumflex over(C)}^(K5),  [eq 2]where {circumflex over (R)}^ is a quantized and scaled version of theimage R; Ĉ^ is a quantized and scaled version of the image C; andparameters K2-K5 are scalar values used to control the relativecontributions of contrast enhancement. K2-K5 are empirically determinedthrough visual inspection.

The above formulation may permit the opaque regions of image (e.g,heart) to be enhanced sufficiently to allow nodules to visualized behindthe heart while simultaneously not “over-enhancing” the air filledregion. Over-enhancement of the air filled region may often lead tonoise enhancement, an undesirable artifact.

The normalized image may be scaled to a target spatial sampling. Therescaled image may be further processed (which may be considered as partof the “enhancement” portion of the noise removal and enhancement) toaccount for localized dark areas introduced (or exaggerated) as part ofthe normalization process. Such further processing may include theaddition of a Laplacian-of-Gaussians (LoG) image to the image. The LoGimage may be clipped and/or scaled so as not to introducediscontinuities.

Various embodiments of the invention may comprise hardware, software,and/or firmware. FIG. 8 shows an exemplary system that may be used toimplement various forms and/or portions of embodiments of the invention.Such a computing system may include one or more processors 82, which maybe coupled to one or more system memories 81. Such system memory 81 mayinclude, for example, RAM, ROM, or other such computer-readable media,and system memory 81 may be used to incorporate, for example, a basicI/O system (BIOS), operating system, instructions/software for executionby processor 82, etc. The system may also include further memory 83,such as additional RAM, ROM, hard disk drives, or othercomputer-readable storage media. Processor 82 may also be coupled to atleast one input/output (I/O) interface 84. I/O interface 84 may includeone or more user interfaces, as well as readers for various types ofstorage media and/or connections to one or more communication networks(e.g., communication interfaces and/or modems), from which, for example,software code may be obtained, e.g., by downloading such software from acomputer over a communication network. Furthermore, other devices/mediamay also be coupled to and/or interact with the system shown in FIG. 8.

The above is a detailed description of particular embodiments of theinvention. It is recognized that departures from the disclosedembodiments may be within the scope of this invention and that obviousmodifications will occur to a person skilled in the art. It is theintent of the applicant that the invention include alternativeimplementations known in the art that perform the same functions asthose disclosed. This specification should not be construed to undulynarrow the full scope of protection to which the invention is entitled.

What is claimed is:
 1. A method of processing an image, comprising: normalizing the image, including resizing the image and adjusting the image bit depth, to obtain a pre-normalized image; obtaining from the pre-normalized image, by an automated processing device, multiple component images, using multi-scale decomposition; performing gray-scale registration based on the multiple component images, wherein the gray-scale registration includes performing histogram matching of the multiple component images in accordance with a target distribution corresponding to at least one histogram obtained based on at least one other image; removing unwanted information from the image using the multiple component images; normalizing the image based on the multiple component images; and enhancing the image by combining a result of the normalizing with a coarse residual component comprising low-frequency information, wherein the combining includes: multiplying the result of the normalizing by a first factor computed based on a quantized and scaled version of the coarse residual component and on a quantized and scaled version of the result of the normalizing; and multiplying the coarse residual component by a scalar design parameter for dynamic range reduction.
 2. The method of claim 1, wherein the multi-scale decomposition comprises at least one decomposition technique selected from the group consisting of: steerable pyramids, Laplacian pyramids, wavelets, and curvelets.
 3. The method of claim 1, wherein the multi-scale decomposition results in an upper-frequency sub-image, and wherein removing unwanted information comprises removing at least one component corresponding to the upper-frequency sub-image.
 4. The method of claim 1, wherein performing a multi-scale decomposition results in a set of multi-scale components that includes the coarse residual component.
 5. The method of claim 4, further comprising performing, on the set of multi-scale components, except for the coarse residual component, at least one of said gray scale registration or said removing unwanted information.
 6. The method of claim 1, further comprising downloading software instructions that, if executed by a processor, cause the processor to perform said normalizing to obtain a pre-normalized image, said obtaining, said performing, said removing, said normalizing based on the multiple component images, and said enhancing.
 7. A non-transitory computer-readable medium containing software instructions that, if executed by a processor, cause the processor to perform operations comprising: normalizing the image, including resizing the image and adjusting the image bit depth, to obtain a pre-normalized image; obtaining from the pre-normalized image, by an automated processing device, multiple component images, using multi-scale decomposition; performing gray-scale registration based on the multiple component images, wherein the gray-scale registration includes performing histogram matching of the multiple component images in accordance with a target distribution corresponding to at least one histogram obtained based on at least one other image; removing unwanted information from the image using the multiple component images; normalizing the image based on the multiple component images; and enhancing the image by combining a result of the normalizing with a coarse residual component comprising low-frequency information, wherein the combining includes: multiplying the result of the normalizing by a factor computed based at least in part on a quantized and scaled version of the coarse residual component and on a quantized and scaled version of the result of the normalizing; and multiplying the coarse residual component by a scalar design parameter for dynamic range reduction.
 8. The medium of claim 7, wherein the multi-scale decomposition comprises at least one decomposition technique selected from the group consisting of: steerable pyramids, Laplacian pyramids, wavelets, and curvelets.
 9. The medium of claim 7, wherein the multi-scale decomposition results in an upper-frequency sub-image, and wherein removing unwanted information comprises removing at least one component corresponding to the upper-frequency sub-image.
 10. The medium of claim 7, wherein performing a multi-scale decomposition results in a set of multi-scale components that includes the coarse residual component.
 11. The medium of claim 10, further comprising performing, on the set of multi-scale components, except for the coarse residual component, at least one of said gray scale registration or said removing unwanted information.
 12. The medium of claim 7, wherein said factor is further computed based on a set of scalar parameters configured to control relative contributions of contrast enhancement.
 13. The method of claim 1, wherein said factor is further computed based on a set of scalar parameters configured to control relative contributions of contrast enhancement. 